clear all
cd "/..."  
use yearly_data.dta, clear
xtset company year

**************************** 定义变量及变量lists *********************************
local all_vars Rev logRev Assets logAst AssetsGrR NumCities LtoAexAd ///
	netDtoE STDitoDi AccRTurnover InvTurnover MABidder SOE SalesExpR FinExpR ///
	AdmExpR AstTurnover ROE ROA 
local w_vars SalesExpR FinExpR AdmExpR AstTurnover OpCashtoRev ROE ROA ///
	AssetsGrR MktShareGr NumCities LtoAexAd netDtoE STDitoDi AccRTurnover InvTurnover
global varlist logRev logRev2 AssetsGrR d.NumCities LtoAexAd STDitoDi ///
	AccRTurnover MABidder SOE i.year
global levarlist l.LtoAexAd lTDTAexAR_sharegro1lvl 
global yvarlist SalesExpR FinExpR AdmExpR AstTurnover ROE ROA
global desvarlist Rev Assets AssetsGrR LtoAexAd STDitoDi AccRTurnover SalesExpR ///
	FinExpR AdmExpR AstTurnover ROE ROA

*按照总资产增速生成虚拟变量
egen assgro3m = median(AssetsGrR_3)
gen assgro3lvl=(AssetsGrR_3>assgro3m) if !missing(AssetsGrR_3)
gen lTDTAexAR_assgro3lvl = l.LtoAexAd * assgro3lvl
gen lAccRTurnover_assgro3lvl = l.AccRTurnover * assgro3lvl

******************************* 描述性统计 ***************************************
*生成相关系数与描述性统计表
logout, save(Tab1_des) word replace: tabstat `all_vars', /// 1.核心变量描述性统计指标
	stats (n mean sd p25 p50 p75) format (%6.2f) column (stats)
	
gen RevSize=(Rev>100) if !missing(Rev)
logout, save(Tab2_byGr) word replace: ttable2 `all_vars', /// 2.大型与小型房企均值差异
	by (RevSize) format (%6.2f) 

logout, save(Tab3_corr) excel replace: pwcorr_a logRev logAst NumCities /// 3.相关系数矩阵
	dNumCities LtoAexAd STDitoDi AccRTurnover MABidder SOE SalesExpR ///
	FinExpR AdmExpR AstTurnover ROE ROA

*绘制核心变量箱型图
graph box logRev if logRev>-2, ///图6，营收对数logRev，-2以下个别极值在图中省略
	over (year, label(labsize(medium) ticks angle(-90) tposition(crossing))) ///
	ytitle(`"{fontface "Times New Roman":营业收入对数 logRev}"', ///
	margin(right) size(large)) marker(1,msize(tiny)) ///
	ylabel(,labsize(medium)) ///
	graphregion(fcolor(white) style(none)) 
	graph display, ysize(2) xsize(4)

graph box LtoAexAd, ///图7，剔除预收款的资产负债率 LtoAexAd
	over (year, label(labsize(medium) ticks angle(-90) tposition(crossing))) ///
	ytitle(`"{fontface "Times New Roman":剔除预收账款的资产负债率 LtoAexAd}"', ///
	margin(right) size(large)) ///
	marker(1,msize(tiny)) ///
	ylabel(,labsize(medium)) ///
	graphregion(fcolor(white) style(none))
	graph display, ysize(2) xsize(4)

graph box AstTurnover, ///图8，总资产周转率 AstTurnover
	over (year, label(labsize(medium) ticks angle(-90) tposition(crossing))) ///
	ytitle(`"{fontface "Times New Roman":总资产周转率 AstTurnover}"', ///
	margin(right) size(large)) marker(1,msize(tiny)) ///
	ylabel(,labsize(medium)) ///
	graphregion(fcolor(white) style(none))
	graph display, ysize(2) xsize(4)

graph box AccRTurnover if AccRTurnover<200, ///图9，应收账款周转率 AccRTurnover，200以上个别极值省略
	over (year, label(labsize(medium) ticks angle(-90) tposition(crossing))) ///
	ytitle(`"{fontface "Times New Roman":应收账款周转率 AccRTurnover}"', ///
	margin(right) size(large)) ///
	marker(1,msize(tiny)) ///
	ylabel(,labsize(medium)) ///
	graphregion(fcolor(white) style(none))
	graph display, ysize(2) xsize(4)

****************************** 基准模型-模型(1)-表3&表4 **********************************
*1a.不加入控制变量的回归-所有样本
foreach yvar in $yvarlist{
	xtreg `yvar' logRev logRev2 AssetsGrR i.year, fe r
	est store reg1a_`yvar'
}
*1b.加入控制变量回归-所有样本
foreach yvar in $yvarlist{
	xtreg `yvar' $varlist, fe r
	est store reg1b_`yvar'
}
*结果汇总-表3
esttab reg1a_SalesExpR reg1b_SalesExpR reg1a_FinExpR reg1b_FinExpR reg1a_AdmExpR ///
	reg1b_AdmExpR reg1a_AstTurnover reg1b_AstTurnover reg1a_ROE reg1b_ROE reg1a_ROA reg1b_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("房企规模效应"),using results_m1.csv,replace

*1c1.分阶段：2001-2007
foreach yvar in $yvarlist{
	xtreg `yvar' $varlist if year<=2007, fe r
	est store reg1c1_`yvar'
}
*结果汇总-表4
esttab reg1c1_SalesExpR  reg1c1_FinExpR  reg1c1_AdmExpR ///
	reg1c1_AstTurnover reg1c1_ROE reg1c1_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("分阶段房企规模效应-2001-2007"),using results_m1c1.csv,replace	

*1c2.分阶段：2009-2015
foreach yvar in $yvarlist{
	xtreg `yvar' $varlist if year>=2008 & year<=2015, fe r
	est store reg1c2_`yvar'
}
*结果汇总-表4
esttab reg1c2_SalesExpR  reg1c2_FinExpR  reg1c2_AdmExpR ///
	reg1c2_AstTurnover reg1c2_ROE reg1c2_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("分阶段房企规模效应-2008-2015"),using results_m1c2.csv,replace	
	
*1c3.分阶段：2016-2021
foreach yvar in $yvarlist{
	xtreg `yvar' $varlist if year>=2016, fe r
	est store reg1c3_`yvar'
}
*结果汇总-表4
esttab reg1c3_SalesExpR  reg1c3_FinExpR  reg1c3_AdmExpR ///
	reg1c3_AstTurnover reg1c3_ROE reg1c3_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("分阶段房企规模效应-2016-2021"),using results_m1c3.csv,replace		
	
*************************** 杠杆率影响-模型（2）-表5 *************************************
*2a.不纳入控制变量的回归
foreach yvar in $yvarlist{
	xtreg `yvar' LtoAexAd l.LtoAexAd lTDTAexAR_assgro3lvl i.year, fe r
	est store reg2a_`yvar'
}
*2b.加入控制变量回归
foreach yvar in $yvarlist{
	xtreg `yvar' l.LtoAexAd lTDTAexAR_assgro3lvl $varlist, fe r
	est store reg2b_`yvar'
}
*结果汇总-表5
esttab reg2a_SalesExpR reg2b_SalesExpR reg2a_FinExpR reg2b_FinExpR reg2a_AdmExpR ///
	reg2b_AdmExpR reg2a_AstTurnover reg2b_AstTurnover reg2a_ROE reg2b_ROE reg2a_ROA reg2b_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("房企杠杆率分析"), ///
	using results_m2.csv,replace

****************************** 周转率影响-模型（3）-表6 **********************************
*3a.不纳入控制变量的回归
foreach yvar in $yvarlist{
	xtreg `yvar' AccRTurnover l.AccRTurnover lAccRTurnover_assgro3lvl i.year, fe r
	est store reg3a_`yvar'
}
*3b.加入控制变量回归
foreach yvar in $yvarlist{
	xtreg `yvar' AccRTurnover l.AccRTurnover lAccRTurnover_assgro3lvl $varlist, fe r
	est store reg3b_`yvar'
}
*结果汇总—表6
esttab reg3a_SalesExpR reg3b_SalesExpR reg3a_FinExpR reg3b_FinExpR ///
	reg3a_AdmExpR reg3b_AdmExpR reg3a_AstTurnover reg3b_AstTurnover ///
	reg3a_ROE reg3b_ROE reg3a_ROA reg3b_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("房企周转率分析"), ///
	using results_m3.csv,replace
	
*********************** 内生性检验-工具变量回归-表7 ****************************
*4a.模型（1）-不纳入控制变量的回归
foreach yvar in $yvarlist{
	xtivreg `yvar' (logRev=L.logRev L2.logRev) i.year, fe first
	est store reg4a_`yvar'
}
*4b.模型（1）-加入控制变量回归
foreach yvar in $yvarlist{
	xtivreg `yvar' (logRev=L.logRev L2.logRev) AssetsGrR LtoAexAd STDitoDi AccRTurnover ///
	MABidder SOE i.year, fe first
	est store reg4b_`yvar'
}
*结果输出汇总-表7
esttab reg4a_SalesExpR reg4b_SalesExpR reg4a_FinExpR reg4b_FinExpR reg4a_AdmExpR ///
	reg4b_AdmExpR reg4a_AstTurnover reg4b_AstTurnover reg4a_ROE reg4b_ROE reg4a_ROA reg4b_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("工具变量回归分析"), ///
	using results_iv.csv,replace
	
****************** 内生性检验-双重机器学习DML-表8（代码见Python文档）*****************
	
*********************** 内生性检验-Oster边界检验-表9 ****************************

*4c.模型（1）-Oster边界检验
foreach yvar in $yvarlist{
	xtreg `yvar' $varlist, fe r level(90)
	global r = e(r2)*1.3 //1.3倍，Oster建议值
	psacalc beta logRev,delta(1) rmax($r) 
	psacalc delta logRev, beta(0) rmax($r)
	est store reg4c_`yvar'
}
*结果输出汇总-表8
esttab reg4c_SalesExpR reg4c_FinExpR reg4c_AdmExpR ///
	reg4c_AstTurnover reg4c_ROE reg4c_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("oster分析"), ///
	using results_oster.csv,replace

************************* 内生性检验-DiD分析-表10*****************************
*按照年份生成政策时间虚拟变量
gen post = (year>2019) //2020"三道红线"政策出台，将数据中2020-2021两个年份划为政策出台后样本
*按照"三条红线"政策生成样本实验组虚拟变量
gen treat_lev = (LtoAexAd>70) if year==2019 
gen treat_netDtoE = (netDtoE>100) if year==2019 
gen treat_cash_STD = (cash_STD>100) if year==2019 
*按照2019年财务数据将红档与橙档房企作为实验组
gen treat_2019 = (treat_lev+treat_netDtoE+treat_cash_STD>=2) if year==2019 
egen treat = max(treat_2019), by (company)
gen treat_post = treat*post //定义did变量
gen treat_post_SOE = treat*post*SOE  //相当于DDD变量，三者交乘项

*国企/民营房企实验组与对照组趋势
**平均债务总额变化趋势
gen debtwithoutAd = totaldebt-Ad
egen DebtnAd_SOEtreat = mean(debtwithoutAd) if SOE==1&treat==1,by (year)
egen DebtnAd_nSOEtreat = mean(debtwithoutAd) if SOE==0&treat==1,by (year)
egen DebtnAd_SOEntreat = mean(debtwithoutAd) if SOE==1&treat==0,by (year)
egen DebtnAd_nSOEntreat = mean(debtwithoutAd) if SOE==0&treat==0,by (year)
//作图
graph set window fontface "Times New Roman"
graph set window fontfacesans "宋体"
graph twoway (scatter DebtnAd_SOEtreat year if company==43,  ///
			mcolor(black) msize(.15cm) connect(1) lcolor(grey) lpattern(solid)) ///
		(scatter DebtnAd_nSOEtreat year if company==53, ///
			msymbol(oh) mcolor(black) msize(.25cm) connect(1) ///
			lcolor(red) lpattern(solid)) ///
		(scatter DebtnAd_SOEntreat year if company==34, ///
			mcolor(black) msize(.15cm) connect(1) lcolor(grey) lpattern(dash)) ///
		(scatter DebtnAd_nSOEntreat year if company==40, ///
			msymbol(oh) mcolor(black) msize(.25cm) connect(1) lcolor(grey) ///
			lpattern(dash)), ///
		legend(label(1 "国企-实验组") label(2 "民企-实验组") ///
			label(3 "国企-对照组") label(4 "民企-对照组")) ///
		legend(size(medium)) ///
		graphregion(fcolor(white) style(none)) ///
		title("样本房企平均负债（剔除预收款）", size(.45cm) color(black)) ///
		xline(2020, lcolor(grey) lwidth(thin) lpattern(dash)) ///
		ttext(2200 2020 "政策发布", size(.35cm)) ///
		ttext(2000 2001.5 "亿元", size(.35cm)) ///
		xlabel(2001(1)2021,angle(vertical)) ///
		xtitle("") ///
		ysize(4) xsize(6) ///
		saving(样本房企平均债务总额变化趋势.png, replace)
	
**平均售额变化趋势
egen Rev_SOEtreat = mean(Rev) if SOE==1&treat==1,by (year)
egen Rev_nSOEtreat = mean(Rev) if SOE==0&treat==1,by (year)
egen Rev_SOEntreat = mean(Rev) if SOE==1&treat==0,by (year)
egen Rev_nSOEntreat = mean(Rev) if SOE==0&treat==0,by (year)
//作图
graph set window fontface "Times New Roman"
graph set window fontfacesans "宋体"
graph twoway (scatter Rev_SOEtreat year if company==43,  ///
			mcolor(black) msize(.15cm) connect(1) lcolor(grey) lpattern(solid)) ///
		(scatter Rev_nSOEtreat year if company==53, ///
			msymbol(oh) mcolor(black) msize(.25cm) connect(1) ///
			lcolor(red) lpattern(solid)) ///
		(scatter Rev_SOEntreat year if company==34, ///
			mcolor(black) msize(.15cm) connect(1) lcolor(grey) lpattern(dash)) ///
		(scatter Rev_nSOEntreat year if company==40, ///
			msymbol(oh) mcolor(black) msize(.25cm) connect(1) lcolor(grey) ///
			lpattern(dash)), ///
		legend(label(1 "国企-实验组") label(2 "民企-实验组") ///
			label(3 "国企-对照组") label(4 "民企-对照组")) ///
		legend(size(medium)) ///
		graphregion(fcolor(white) style(none)) ///
		title("样本房企平均销售额合计",size(.45cm) color(black)) ///
		xline(2020, lcolor(grey) lwidth(thin) lpattern(dash)) ///		
		xtitle("") ///
		xlabel(2001(1)2021,angle(vertical)) ///
		ttext(670 2020 "政策发布", size(.35cm)) ///
		ttext(600 2001.5 "亿元", size(.35cm)) ///
		ysize(4) xsize(6) ///
		saving(样本房企平均销售总额变化趋势.png,replace)
		
**平均总利润额变化趋势
egen totprofit_SOEtreat = mean(totalprofit) if SOE==1&treat==1,by (year)
egen totprofit_nSOEtreat = mean(totalprofit) if SOE==0&treat==1&company!=18,by (year)
	*计算利润总额时剔除了异常值中国恒大，其2021年度亏损额过高
egen totprofit_SOEntreat = mean(totalprofit) if SOE==1&treat==0,by (year)
egen totprofit_nSOEntreat = mean(totalprofit) if SOE==0&treat==0,by (year)
//作图
graph set window fontface "Times New Roman"
graph set window fontfacesans "宋体"
graph twoway (scatter totprofit_SOEtreat year if company==43,  ///
			mcolor(black) msize(.15cm) connect(1) lcolor(grey) lpattern(solid)) ///
		(scatter totprofit_nSOEtreat year if company==53, ///
			msymbol(oh) mcolor(black) msize(.25cm) connect(1) ///
			lcolor(red) lpattern(solid)) ///
		(scatter totprofit_SOEntreat year if company==34, ///
			mcolor(black) msize(.15cm) connect(1) lcolor(grey) lpattern(dash)) ///
		(scatter totprofit_nSOEntreat year if company==40, ///
			msymbol(oh) mcolor(black) msize(.25cm) connect(1) lcolor(grey) ///
			lpattern(dash)), ///
		legend(label(1 "国企-实验组") label(2 "民企-实验组") ///
			label(3 "国企-对照组") label(4 "民企-对照组")) ///
		legend(size(medium)) ///
		graphregion(fcolor(white) style(none)) ///
		title("样本房企平均总利润合计",size(.45cm) color(black)) ///
		xline(2020, lcolor(grey) lwidth(thin) lpattern(dash)) ///		
		xtitle("") ///
		xlabel(2001(1)2021,angle(vertical)) ///
		ttext(109 2020 "政策发布", size(.35cm)) ///
		ttext(103 2001.5 "亿元", size(.35cm)) ///
		ysize(4) xsize(6) ///
		saving(样本房企平均利润总额变化趋势.png,replace)
		
*DiD回归
xtreg logRev treat_post treat_post_SOE AssetsGrR LtoAexAd STDitoDi AccRTurnover ///
	MABidder i.year, fe r
	est store reg5a_logRev
xtreg AssetsGrR treat_post treat_post_SOE logRev LtoAexAd STDitoDi AccRTurnover ///
	MABidder i.year, fe r
	est store reg5a_AssetsGr
foreach yvar in $yvarlist{
	xtreg `yvar' treat_post treat_post_SOE logRev AssetsGrR LtoAexAd STDitoDi ///
	AccRTurnover MABidder i.year, fe r
	est store reg5a_`yvar'
}
**结果汇总
esttab reg5a_logRev reg5a_AssetsGr reg5a_SalesExpR reg5a_FinExpR reg5a_AdmExpR ///
	reg5a_AstTurnover reg5a_ROE reg5a_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("DiD回归"), ///
	using results_5a.csv,replace

********************* 稳健性检验-替换核心解释变量-***********************************
*准备工作：生成新的变量list
global varlist_rc logAst logAst2 RevGrR d.NumCities netDtoE STDitoDi ///
	InvTurnover MABidder SOE i.year

*6a.稳健性检验-模型（1）替换解释变量
foreach yvar in $yvarlist{
	xtreg `yvar' logAst logAst2 RevGr i.year, fe r
	est store reg6a_`yvar'
}
*6b.稳健性检验-模型（1）替换解释变量-加入控制变量
foreach yvar in $yvarlist{
	xtreg `yvar' $varlist_rc, fe r
	est store reg6b_`yvar'
}
*结果汇总
esttab reg6a_SalesExpR reg6b_SalesExpR reg6a_FinExpR reg6b_FinExpR reg6a_AdmExpR ///
	reg6b_AdmExpR reg6a_AstTurnover reg6b_AstTurnover reg6a_ROE reg6b_ROE reg6a_ROA ///
	reg6b_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("房企规模效应_稳健性检验_指标替换"), ///
	using results_rc1_m1.csv,replace

*7a.稳健性检验-模型（2）替换解释变量
gen lnetDtoE_assgro3lvl = l.netDtoE * assgro3lvl 
foreach yvar in $yvarlist{
	xtreg `yvar' netDtoE l.netDtoE lnetDtoE_assgro3lvl i.year, fe r
	est store reg7a_`yvar'
}
*7b.稳健性检验-模型（2）替换解释变量-加入控制变量
foreach yvar in $yvarlist{
	xtreg `yvar' l.netDtoE lnetDtoE_assgro3lvl $varlist_rc, fe r
	est store reg7b_`yvar'
}
*结果汇总
esttab reg7a_SalesExpR reg7b_SalesExpR reg7a_FinExpR reg7b_FinExpR reg7a_AdmExpR ///
	reg7b_AdmExpR reg7a_AstTurnover reg7b_AstTurnover reg7a_ROE reg7b_ROE ///
	reg7a_ROA reg7b_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("房企杠杆率分析_稳健性检验_指标替换"), ///
	using results_rc1_m2.csv,replace

*8a.稳健性检验-模型（3）替换解释变量
gen lInvTurnover_assgro3lvl = l.InvTurnover * assgro3lvl
foreach yvar in $yvarlist{
	xtreg `yvar' InvTurnover l.InvTurnover lInvTurnover_assgro3lvl i.year, fe r
	est store reg8a_`yvar'
}
*8b.稳健性检验-模型（3）替换解释变量-加入控制变量
foreach yvar in $yvarlist{
	xtreg `yvar' l.InvTurnover lInvTurnover_assgro3lvl $varlist_rc, fe r
	est store reg8b_`yvar'
}
*结果汇总
esttab reg8a_SalesExpR reg8b_SalesExpR reg8a_FinExpR reg8b_FinExpR reg8a_AdmExpR ///
	reg8b_AdmExpR reg8a_AstTurnover reg8b_AstTurnover reg8a_ROE reg8b_ROE ///
	reg8a_ROA reg8b_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("房企周转率分析_稳健性检验_指标替换"), ///
	using results_rc1_m3.csv,replace

************************ 稳健性检验-剔除小规模样本-**********************************

*分组
gen above = 0
egen median_logRev = median(logRev)
replace above = 1 if logRev >= median_logRev

*9a.稳健性检验-模型（1）剔除小规模样本
foreach yvar in $yvarlist{
	xtreg `yvar' logRev logRev2 AssetsGrR i.year if above==1, fe r
	est store reg15_`yvar'
}
*9b.稳健性检验-模型（1）剔除小规模样本-加入控制变量
foreach yvar in $yvarlist{
	xtreg `yvar' $varlist if above==1, fe r
	est store reg16_`yvar'
}
*结果汇总
esttab reg15_SalesExpR reg16_SalesExpR reg15_FinExpR reg16_FinExpR reg15_AdmExpR ///
	reg16_AdmExpR reg15_AstTurnover reg16_AstTurnover reg15_ROE reg16_ROE ///
	reg15_ROA reg16_ROA, ///
	b(%9.3f) star(* 0.1 ** 0.05 *** 0.01) se scalars(N r2 F p) mtitles ///
	title("房企规模效应_稳健性检验_样本缩减"), ///
	using results_rc2_m1.csv,replace

************************ 稳健性检验-置换检验-**********************************
permute treat_post beta=_b[treat_post], reps(1000): xtreg logRev treat_post ///
	treat_post_SOE AssetsGrR LtoAexAd STDitoDi AccRTurnover ///
	MABidder i.year, fe r
	est store reg5a_logRev
permute treat_post beta=_b[treat_post], reps(1000):xtreg AssetsGrR treat_post ///
	treat_post_SOE logRev LtoAexAd STDitoDi AccRTurnover ///
	MABidder i.year, fe r
	est store reg5a_AssetsGr
foreach yvar in $yvarlist{
	permute treat_post beta=_b[treat_post], reps(1000):xtreg `yvar' treat_post ///
	treat_post_SOE logRev AssetsGrR LtoAexAd STDitoDi ///
	AccRTurnover MABidder i.year, fe r
}




